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RESUMEN 

Se estudian los estados finales dc cquilibrio que surgen de la inestabilidad de Parker cuando partimos de una 
configuration cih'ndrica de gas en cquilibrio magnctohidrostatico en un campo gravitational radial y con un 
campo magnctico longitudinal. Nucstro objetivo es comparar los estados de equilibrio no lincales con los que 
se obtienen en los sistemas con geometrfa Cartesiana. Se presentan los mapas de densidad y dc las lincas de 
campo magnctico en ambas geomctrias para un campo gravitational dc intensidad constantc. Encontramos 
que el fiotamiento magnctico es menos eficicnte bajo simetria axial que en una atmosfera Cartesiana. Como 
consecucncia, las condensaciones que se forman en el modelo axisimetrico tienen menor densidad columnar. 
Por ende, el cociente entre la presion magnetica y termica en el estado final toma valores mas extremos bajo 
simetria Cartesiana. Se discuten tambien algunos modelos en los que el campo gravitational no es uniforme. 

ABSTRACT 

We study the final equilibrium states of the Parker instability arising from an initially unstable cylindrical 
equilibrium configuration of gas in the presence of a radial gravitational field and a longitudinal magnetic field. 
The aim of this work is to compare the properties of the nonlinear final equilibria with those found in a system 
with Cartesian geometry. Maps of the density and magnetic field lines, when the strength of the gravitational 
field is constant, are given in both geometries. Wc find that the magnetic buoyancy and the drainage of gas 
along field lines are less efficient under axial symmetry than in a Cartesian atmosphere. As a consequence, 
the column density enhancement arising in gas condensations in the axially-symmetric model is smaller than 
in Cartesian geometry. The magnetic-to-gas pressure ratio in the final state takes more extreme values in the 
Cartesian model. Models with non-uniform radial gravity are also discussed. 
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1. INTRODUCTION 

Parker (1966) demonstrated that a layer that is 
supported against the vertical gravity by the pres- 
sure of interstellar gas, horizontal magnetic field and 
cosmic rays, is always unstable to long wavelength 
deformations of the field lines because of the buoy- 
ancy of magnetic field and cosmic rays. Much work 
has been done to understand the role of the Parker 
instability in modeling the structure of the interstel- 
lar medium (ISM) and, more specifically, its role in 
gathering gas to trigger the formation of giant molec- 
ular clouds (e.g., Mouschovias et al. 1974; Blitz & 
Shu 1980; Mouschovias ct al. 2009). Giz & Shu 
(1993), Kim ct al. (1997) and Kim & Hong (1998) ex- 
amined how the non-uniform nature of the Galactic 
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gravity might affect the length and timescales of the 
Parker instability. Many studies have been devoted 
to investigate the effect of other physical ingredients 
on the Parker instability (Shu 1974; Nakamura et 
al. 1991; Hanawa et al. 1992; Hanasz & Lesch 2000 
Kim ct al. 2000; Santillan ct al. 2000; Kim et al. 2001 
Franco ct al. 2002; Kosihski & Hanasz 2006, 2007 
Lee & Hong 2007; Mouschovias ct al. 2009). 

The structure of the final static state result- 
ing from the Parker instability was derived by 
Mouschovias (1974, hereafter M74) in the Carte- 
sian model, namely, for a non-rotating layer under 
a vertical gravity, with a horizontal magnetic field 
in the initial equilibrium state. In the unperturbed 
state, the gas is in hydrostatic equilibrium in a strat- 
ified layer, so that all the quantities depend only on 
the distance from the midplane y, and the magnetic 
field is plane-parallel and runs along the ^-direction. 
In this two-dimensional Cartesian problem, one as- 
sumes that all the quantities are independent of the 
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third dimension, z, at any time. In the final state, 
the magnetic field rises in certain regions and sink in 
others forming bulges and valleys. Matter loads onto 
the field lines by draining down from the top region 
of the bulge and sinking into the valleys. So far, the 
final states of the Parker instability have been dis- 
cussed in Cartesian models. In various astrophysical 
systems, however, the gas tends to form structures 
with other symmetries, such as filaments or elon- 
gated structures (e.g., Genzel & Stutzki 1989; Alfaro 
ct al. 1992; Ryu et al. 1998; Conselice et al. 2001; 
Salome et al. 2006; Jackson et al. 2010). Suppose 
now that we have a very elongated self-gravitating 
filamentary structure with a magnetic field travel- 
ing along the major axis of this structure, such that 
it has the effective geometry of a infinite cylinder. 
If the filament is supported against the radial self- 
gravity of the system by the pressure of gas and the 
longitudinal magnetic field, long wavelengths defor- 
mations of the magnetic field lines will allow matter 
to slide radially along field lines under the action of 
the gravitational field towards the symmetry axis, 
whereas the magnetic field tends to expand in ra- 
dius, through a process of magnetic buckling (see 
Fig. 1). If the disturbances are axisymmetric, we 
expect the alternate formation of magnetic bottle- 
necks in the condensations of gas separated by por- 
tions of inflated magnetic lines or arcades, resem- 
bling a sausage mode. The final state is expected 
to be different than in a vertically stratified disk as 
those models currently used to study the Parker in- 
stability in galactic disks. In fact, the response of 
the system and the final equilibrium state depend 
on the adopted model: Cartesian or axisymmetric. 

In this paper, we study the nonlinear final equi- 
libria of the Parker instability of a non-rotating, 
cylindrical initial state. In the initial state, the mag- 
netic field is assumed to be parallel to the axis of 
the cylinder, which is taken to coincide with the z- 
axis. The cylinder is supported against the action 
of a radial gravitation field gn by the thermal and 
magnetic pressures (see Fig. 1). We will assume that 
the gas is isothermal and evolves under flux-freezing 
conditions. For axisymmetric perturbations, all the 
quantities only depend on z and R (not on the az- 
imuthal angle (j>) and hence the problem is also two- 
dimensional. The evolution of the Parker instability 
in the (z, _R)-plane is expected to be qualitative anal- 
ogous in the main features, to the evolution found in 
the (x, y)-plane for the case of a plane-parallel ini- 
tial state by M74. Throughout the paper, we will use 
the coordinate z to denote the distance on the axis 
of symmetry in the axisymmetric model and should 




Fig. 1. Sketch of the magnetic field configuration in 
the initial (left panel) and final (right panel) equilibrium 
states in axial symmetry. 

not be confused with the third component of the 
Cartesian model. 

The aim of this work is to compare in detail how 
the final distributions of density and magnetic field 
vary from Cartesian to axisymmetric geometry 3 . Ac- 
curate predictions of the nonlinear final states, an 
interesting problem in itself, can be very useful for 
the purpose of testing MHD codes. We will follow 
Mouschovias' procedure (1974, 1976) to construct 
the final equilibrium states. Nakano (1979) and 
Tomisaka et al. (1988) used the same technique to 
find the structure of axisymmetric magnetized clouds 
that can be reached from nonequilibrium spherical 
configurations. By contrast, we start from a cylin- 
drical configuration in magnetohydrostatic equilib- 
rium as those used to study the gravitational col- 
lapse and fragmentation of filamentary clouds (e.g., 
Nakamura et al. 1995). In these simulations, how- 
ever, the gravitational instabilities play a dominant 
role and not the Parker instability itself. In order 
to isolate the outcomes of the Parker instability, we 
will not consider self-gravity of the gas. In addi- 
tion, we will also ignore the destabilizing influence 
of the interchange mode and the magnetic field cur- 

3 Note that the change in the geometry of the system is, by 
no means, equivalent to a coordinate transformation. 
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vature in the equilibrium states (Asseo et al. 1978, 
1980; Lachieze-Rey et al. 1980), which could lead to 
the development of substructure. The study of these 
instabilities requires a full magnetohydrodynamical 
model and they will not be discussed further here. 

Our paper is organized as follows. In Section 
2, we describe the initial equilibrium configurations. 
The elliptic differential equations that give the equi- 
librium states under flux-freezing conditions are pro- 
vided in Section 3. In Section 4 we compare the 
resulting features of the final states arising from 
the Parker instability in Cartesian and axisymmetric 
models. A summary of the results are presented in 
Section 5. 

2. INITIAL STATES 

We consider first the initial state in the Carte- 
sian model, that is, a planar layer that is supported 
against the vertical gravity field by the gas pressure 
P and a horizontal magnetic field. We will use a 
Cartesian coordinate system (x,y), where x is the 
horizontal coordinate and y is the distance to the 
midplane of the layer, so that B = (B 0yX (y),0) in 
the initial state. Throughout the paper, the sub- 
script signifies the initial state. Magnetohydro- 
static equilibrium is satisfied if 

V (p + = -pV1>, (1) 

where ip = tp(y) is the gravitational potential. If 
the gas is isothermal with thermal sound speed c s 
and the ratio of the magnetic to gas pressures, a = 
Bq/SitPq, is constant, the equilibrium equation reads 

(l + a^Vlnpo = -W, (2) 
which has the solution 

Po(y) = Po(0)exp f- ( 1+ ^ a ) c 2 j ' ( 3 ) 

with ip(0) = 0. If the gravitational field g = —S7ip is 
taken to be vertical and constant (but it reverses its 
direction across the midplane), we have ip = g\y\, 
with g = constant > 0. In this particular case, 
the mass density decays exponentially with a typ- 
ical scale H of (1 + a)c^/g, whereas the magnetic 
field strength decays with a characteristic scale 2H. 

If, instead of a planar atmosphere, we assume 
that the gas extends to infinity along a cylinder 
whose axis coincides with the z-axis, and is threaded 
by an axial magnetic field B z , the density profile at 



equilibrium is: 

Mi?)=Po(0)cxp(- ir ^- f ), (4) 

where R is the radial distance in cylindrical coordi- 
nates and ip(R) is the gravitational potential, with 
V'(O) = 0- In fact, wc assume that the gravitational 
acceleration has only a radial component. If the ra- 
dial gravitational force is taken to be constant, it 
holds that ip = gR. Therefore, the radial scalelength 
for the density is H = (1 + a)c 2 s / '<?, as in Cartesian 
geometry. 

We see that if the gravitational potential ip{y) in 
the infinite planar layer is identical to the gravita- 
tional potential ip(R) in the cylindrical model, the 
profiles of density, magnetic field and gas pressure 
along y and R will be identical in both models. How- 
ever, these systems are no longer equivalent when the 
Parker instability arises. In the Cartesian model, 
magnetic field lines are allowed to inflate only in the 
vertical direction, whereas the expansion of field lines 
as well as the contraction of matter occur radially in 
the axisymmetric model (see Fig. 1). The question 
that arises is how the distributions of density, pres- 
sure and field lines in an azimuthal section [i.e. in the 
(z, i?)-plane] look like as compared to those found in 
a vertical section [i.e. in the (x, y)-planc] in the well- 
studied Cartesian model. 

When studying the Parker instability at the so- 
lar neighborhood, it is common to use a local Carte- 
sian frame where the azimuthal direction is identified 
with the horizontal coordinate x and the vertical di- 
rection of the Galaxy with y. In the two-dimensional 
approximation, the third radial direction is ignored. 
As a cautionary comment, our coordinate R used 
throughout the paper should not be confused with 
the radial direction of the Galaxy. 

3. THE EQUILIBRIUM EQUATIONS 
If a disturbance of sufficiently long wavelength is 
applied along the initial magnetic field to the initial 
equilibrium states described in the previous Section, 
the Parker instability develops in time and parts of 
the magnetic field lines bulge outward until the ten- 
sion of the field lines increases sufficiently for an equi- 
librium state to be established. In the final states 
arising from the Parker instability, the magnetic field 
lines are curved and the magnetohydrostatic force 
equation is 

-VP- pVtp + -Jx B = 0. (5) 

c 

where J is the current density. The coordinate sym- 
metries of the unperturbed states (see Section 2) re- 
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duce the problem to the (x, y)-plane in the Cartesian 
model and to the (z, iZ)-plane in the axisymmetric 
model. Under these symmetries, flux-freezing allows 
to transform the magnetostatic equations of equi- 
libria into a second-order elliptic partial differential 
equation (Dungcy 1953), conventionally called the 
Grad-Shafranov equilibrium equation. The beauty 
of this approach is that the final equilibrium states 
can be found with no need to solve a time-dependent 
problem (M74; Mouschovias 1976). In the following 
we give the Grad-Shafranov equilibrium equations 
and the boundary conditions. Details on its deriva- 
tion can be found in M74, Mouschovias (1976) and 
Nakano (1979). 

3.1. The Cartesian model 
In the Cartesian layer with one ignorable coor- 
dinate, so that all quantities depend only on (x,y), 
the magnetic vector potential can be written as A = 
A(x, y)e z and then B x = dA/dy and B y = —dA/dx. 
One can show that the scalar function q(x,y), de- 
fined by 



q = Pexp 



(0) 



is constant on a field line at magnetohydrostatic equi- 
librium, and is given by: 



q{A) = 



2 dA 



x * /2 , dy(x,A) 
dx — — — exp 
a A 



%l>{x,A) 



(7) 

where X x is the perturbation wavelength along the 
initial magnetic field and dm is the mass per unit 
length (along the ignorable coordinate) in a flux tube 
between field lines characterized by A and A + dA 
and between x = and x = X x (i.e. the mass-to-flux 
ratio). Under flux- freezing conditions, the mass-to- 
flux ratio is a constant of motion and, therefore, can 
be determined from the initial equilibrium configu- 
ration: 



dm 
~dA 



= A. 



Po(A) 
B (AY 



(8) 



Finally, the magnetic equilibrium equation (5) 
may be written in terms of A as 



V Z A 



A d 1 

-4:TT—— exp 

dA \ c 



(9) 



(Dungcy 1953; M74). Once the boundary conditions 
are specified, Equations (7) and (9) can be solved 
numerically by using a iterative scheme (M74). 

The boundary conditions we take are the same as 
those in M74. The system is assumed to be periodic 



in x and symmetric about the x-axis so that Equa- 
tion (9) is solved in the rectangle < x < X max with 
AT max — X x /2 and < y < Y max , with the boundary 
conditions 



dA 

dx 



0. 



(10) 



=o,x n 



The upper boundary must be far enough from the 
midplane not to affect the evolution since a small 
vertical size of the computational domain suppresses 
the instability (e.g., M74; Mouschovias 1996). At the 
upper boundary we impose that the field lines are 
not deformed, A(x,Y m£lx ) = A a (Y max ). Due to the 
imposed reflection symmetry about the midplane, 
the field line originally coinciding with the a;-axis re- 
mains undeformed. 

In the following, we provide the equations for the 
particular case of a constant external field. Although 
a uniform acceleration is only realistic at high alti- 
tudes from the midplane, we will concentrate on this 
case for simplicity and to facilitate comparison with 
some previous work that use a constant gravitational 
field (e.g., M74; Mouschovias et al. 2009). A realis- 
tic model of the gravitational potential is important 
to derive the growth rate of the instability and the 
most unstable wavelength, which determines the sep- 
aration of the magnetic valleys (Kim & Hong 1998; 
Section 4.2). For now, we are mainly interested in 
the main conceptual features that arise in the final 
equilibria state when the geometry of the problem 
is changed, even if the adopted spatial dependence 
of the gravitational field lacks a strong astrophysical 
motivation. 

As said in Section 2, if the external gravita- 
tional field is constant, the gravitational potential 
is ip = g\y\, with g a positive and constant acceler- 
ation. The equations are put in dimensionless form 
by choosing the following natural units: c 2 s / g, c s and 
po(0), for length, velocity and density, respectively, 
The magnetic vector potential will be measured in 
units of —2HBq(Q). Note that the dimensionless 
vertical scale is H = (1 + a) in these units. The 
dimensionless of Equation (9) has the form 



V 2 A = Q xy (x,y;a), 



(11) 



where 



dq 



^ xy 8a(l + a) 2 dA 
In these units, 



exp(-j/). (12) 



... 1 dm 
qiA) = 2dA 



(i 



A - /2 , dy(x,A) , . ... 
dx yy dA ' exp[-y(x,A)} 



(13) 
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where 



dm 
~dA 



-2X X (1 + a) A. 



(14) 



Under constant gravity condition, the vector po- 
tential in the plane-parallel unperturbed state is 



A (y) = cxp 



V 



2(1 + a) 

and the source term of Eq. (11) is 
1 ( V 



Q 



(i) 



4(1 + a y 



■ cxp 



2(1 + a) 



(15) 



(16) 



In this initial state, the relation between A and q is 
A (q) = q- 1/{2a) (see also M74). 

3.2. The axisymmetric model 
Now, we consider a three-dimensional geometry 
with axial symmetry about the z-axis and longitu- 
dinal magnetic field (Bf = 0). The magnetic field 
is related to the potential vector A — A(z, R)e$ 
by B z = R- 1 d(RA)/dR and B R = -dA/dz. Un- 
der axisymmetric disturbances, the magnetic flux 
$(i?, z) = RA is constant on a magnetic surface and 
is also a constant of the motion. In this geometry, 
the function q defined in Eq. (7) is given by 



Az/2 



dR 2 (z,^>) 



exp 



(17) 

where X z is the longitudinal perturbation wavelength 
and dm is the mass in a flux tube between $ and 
$ + d$ and between z — to z — X z . Again, q is 
a constant along field lines in equilibrium configura- 
tions. Under flux- freezing conditions, the mass-to- 
flux ratio is a constant of the motion, which can be 
determined from the initial configuration by: 



dm 



-DO 



(18) 



The balance of forces (Eq. 5) in the final magnetic 
configuration can be rewritten in terms of $ as 



2 <9$ 
RdR 



A D 2 d 1 

-vkR — cxp 
cf<& 



(19) 



(Mouschovias 1976; Nakano 1979). Equation (19) to- 
gether with Eq. (17) enable us to derive $(z, R) once 
the boundary conditions are specified. The mag- 
netic field can be found just by simple derivation: 
B z = R- l d$/dR and B R = -R- x d$/dz. 

The computational domain is a rectangular box 
of < z < Z max and < R < i? max - Because of the 



symmetries of the initial and final states, the length 
of the computation box in the ^-direction, Z max , is 
taken to be A z /2, half the wavelength of the pertur- 
bation. Periodicity in z is translated to 



<9$ 

dz 



0. 



(20) 



z=0,Z n 



The definition of $ means that $(0, z) = at any 
time. Assuming that the deformation of field lines 
can be neglected on the surface of the cylinder, we 
impose that $(z,i? max 

The equivalent condition of constant gravity is to 
assume that the radial gravitational field is constant, 
i.e. g = —ge R , with g a positive constant accelera- 
tion. We may use the same units for distance, ve- 
locity and density as those described in the previous 
subsection. The magnetic flux will be measured in 
units of — 4H 2 Bq(Q). Equation (19) in dimcnsionless 
form is: 



2 <9$ 
RdR 



Qzr(z,R; a), 



(21) 



with 



and 



izR 



1 dm 
o($) = 



R 2 dq 
32a(l + a) 4 d¥ 

x ' /2 dR 2 {z,<S>) 

dz 

o 5$ 



exp(-ii), (22) 



exp[-i?(z,$)] 



where the dimcnsionless mass-to-flux is 
dm 



?7r(l + a) 2 \ z exp 



2(1 + a) 



(23) 

(24) 
(25) 



By comparing Eqs. (21)-(22) with Eqs. (11)-(12), we 
see that the dynamical variable $ also satisfies a non- 
linear elliptical differential equation, but with a dif- 
ferent differential operator and source term than the 
equation governing A in the Cartesian model. 

In order to gain some physical insight, it is useful 
to compute the function in the initial equilib- 
rium configurations which, of course, are solutions 
of the differential equation (21). In a model with 
uniform radial gravitational field, the magnetic flux 
and the g-function in the unperturbed state are, re- 
spectively: 



MR) 



i 



R 



2(1 + a) 



cxp 



R 



2(1 + a). 



(26) 
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Fig. 2. Azimuthal section of the magnetic flux &(z, R) in 
the axially-symmetric final state for a = 1 and gn = g = 
const, when the upper boundary is placed at i? max = 24 
(top) and at i? max = 33.5 (middle). Magnetic field lines 
coincide with the isocontours of <E>. The color bar is in 
natural logarithmic scale. The unit length is <? s /g. A 
comparison of the radial profile of the magnetic flux along 
cuts at z — and z = Z max is given in the lower panel 
when the boundary is at 7? max = 24 and i? max = 33.5. 




log(A-') log(r') 



Fig. 3. The shape of the function q(A) (left panel) and 
g($) (right panel) in the initial state (dashed lines) and 
in the final state (solid lines), with a = 1 and constant 
gravity, for X max = 15 and Z max = 15, respectively. Note 
that the range of the g-axis is different. 

and 

q (R) = cxp (jf^ ■ (27) 

At small enough R, say R -c 2(1 + a), they can be 
approximated by $o(R) ~ exp(— i?/2[l + a}) and 
<7o($o) = ^o^ 2 "- Interestingly, they have exactly 
the same forms as Ao(y) and qo (A) in the Cartesian 
model, respectively (see Eq. 15 and the very end of 
Section 3.1). At large radii, R ^> 2(1 + a), or equiv- 
alcntly ]nq 3> 2a, the relationship between $o a- n d 
qo is $o — (2a)~ 1 9o 1 In go- Therefore, at these 
large q values, q increases with $ _1 slightly faster 
than the rate at which q increases with A -1 in the 
Cartesian model (q = A~ 2a ). In Section 4, we will 
calculate the dependence of the function q on $ in 
the final axially-symmetric equilibrium state. The 
shape of <?($) as compared to q(A) in the Cartesian 
model will be also discussed. 

The source term of the differential equation in 
the initial cylindrical configuration with constant gn 
is: 

Ql ^w^ R ^(-^j)- (28) 

By comparing Eqs. (16) and (28), we see that Q^ R 
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Fig. 4. The dependence of the functions Q xy and Q z r on y and R, respectively in the initial states (dashed lines) and 
in the final states (solid lines), with a = 1 and constant gravity. In the left panels, Q xy is shown along the lines x = 
and x = X max for two different wavelengths: X x /2 = 9 and 15. The right panels display Q z r along the lines z = and 
z — Zmax for A z /2 = 9 and 15, as quoted at the corners of these panels. 
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has a similar dependence on R as Qxy on y, except 
by the factor R/(2 + 2a). As a consequence, Q^ R 

has a maximum at 2(1 + a) = 2H whereas Qxy is a 
monotonically decreasing function of y. In order to 
facilitate comparison between both cases, it is con- 
venient to use Q xy defined as Q xy = yQ xy /(2 + 2a) 
because both Q xy and Q z r are identical at the ini- 
tial state. However, they are expected to depart one 
from each other in the nonlinear equilibrium state 
because of the changes introduced by the different 
geometries. The form of Q xy and Q z r in the final 
equilibrium states will be discussed in the next Sec- 
tion. 



4. FINAL EQUILIBRIUM STATES 

The final equilibrium state can be reached from 
the initial configuration through continuous defor- 
mations of field lines using a iterative method as de- 
scribed in M74. To do so, we add a perturbation 
to the magnetic variables in the initial equilibrium 
state, 



SA(x, y) = -fiA (y) sin ( ) cos 



(29) 



for the Cartesian model and 



6$(z,R) = -/z$ (-R)sin 



ttR 



z„ 



(30) 

for the axisymmetrical model. Here [i is the am- 
plitude of the perturbation which we chose between 
0.02 and 0.06. Note that the perturbations are sym- 
metric about y = in the Cartesian model 4 . 

Our main interest is to compare the final equilib- 
rium states in Cartesian and axisymmetrical models 
using the same external gravitational potential along 
y and R, respectively. In this way, the comparison 
can be done on a common ground. In the next sub- 
section, wc will do so for models with constant grav- 
itational acceleration. A more realistic gravitational 
potential will be studied in Section 4.2 in order to ex- 
amine how it might affect the final nonlinear states. 

4.1. Constant gravitational acceleration 

Parker (1966) showed that a plane-parallel atmo- 
sphere with uniform gravity and a horizontal mag- 
netic field is unstable to the undular mode if the 
horizontal wavelength of the perturbation is larger 
than a critical value, namely, X x > A cr it = 47r(l + 



4 The Cartesian model allows also solutions where field lines 
cross the midplane. They are called odd-parity solutions or 
midplane antisymmetric modes. 



a)(2a + 1)~ 1//2 . For the reference value a = 1, the 
dimensionlcss critical wavelength is 14.5 or, cquiv- 
alently, 7.26H. Magnetic tension can stabilize dis- 
turbances with shorter wavelengths. The horizontal 
wavelength of the mode with the maximum growth 
rate depends on the vertical wavelength. For a typ- 
ical value of i^nax = 24, the maximum growth rate 
occurs at X max = 12 (i.e. X x = 24). 

The stability of magnetized cylinders has been 
considered in the past (e.g., Stodolkicwicz 1963; Na- 
gasawa 1987) but, to our knowledge, none derived 
the dispersion relation under the assumption of con- 
stant gravitational field. Nevertheless, the linearized 
equations of motion under axial geometry and ax- 
isymmctric perturbations are identical to eqs. (III- 
3)- (III-6) quoted in Parker (1966), once the sub- 
scripts y and z in Parker's notation are replaced for 
z and R, except for the curvature term vr/R that 
appears as an extra term in the continuity and en- 
ergy equations. Since the curvature term becomes 
small at scales comparable or larger than A cr it, the 
critical wavelength is essentially the same in both 
geometries. In fact, we have checked that our itera- 
tive scheme converges to the initial state as long as 
X z < 7.5-ff (when a = 1 is used), indicating that 
they are stable. 

In the Cartesian atmosphere, M74 found that if 
the upper boundary is far enough from the galactic 
plane, i.e. Y ms ,x 3> (l + a ) m dimensionless units, the 
evolution of most of the matter in the domain is not 
affected by the location of the upper boundary. The 
reason is that more than 90% of the energy of the 
initial state resides under y < 3.5(1 + a) (sec M74). 
In the axisymmetrical model, we can apply the same 
reasoning to expect that the solution is not signifi- 
cantly altered as long as the box size is much larger 
than the effective scalchcight, defined as the radius R 
at which <3?o decays by a factor e, namely, 4.3(1 + a) 
[note that the effective scaleheight of Aq in the Carte- 
sian model is 2(1 + a)]. In Figure 2 we compare the 
the magnetic flux of the final state in the axisymmet- 
rical model using a = 1 and A 2 /2 = 12, when the up- 
per boundary, i? max , is located at 24 and at 33.5. We 
see that both maps are almost identical at R < 20. 
The main difference is that upper magnetic field lines 
in the magnetic bulge region inflate further out when 
the boundary is far away because a larger column of 
material is unloaded and brought down to the valley 
of the curved magnetic field lines. Therefore, this 
upper region becomes lighter and more buoyant. By 
contrast, cuts of <& along the valley z = are very 
similar in both cases (see the lower panel of Fig. 2). 
This means that the extra material that has been 
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loaded into the valley when the upper boundary is 
placed at 33.5, does not cause significant additional 
compression in the radial direction. 

In order to have the same dynamical range in 
all our calculations, we will place our upper bound- 
ary at the distance at which the variable of interest 
(A for the Cartesian model and $ for the axisym- 
metric case) decreases to exp(— 6.12) = 0.0022 of 
its maximum value. With this convention, Y max = 
12.2(1 + a) and i? max = 16.7(1 + a). Interestingly, 
in the axisymmetric case, the corresponding density 
at i? max is only 5.4 x 10 -8 of its value at R = 0. 

Different magnetostatic equilibrium configura- 
tions can be found by choosing the function q{A) for 
the Cartesian model and q(&) for the axisymmetric 
one. Parker (1966) explored the cases in which q{A) 
is cither a linear or quadratic function of A. In prac- 
tice, these functions are not free because they are 
defined by the initial conditions. M74 emphasized 
that q increases when A decreases in the Cartesian 
model. In Fig. 3, the functions q(A) for the Cartesian 
layer and q{&) for the axially-symmetric configura- 
tion are plotted for a = 1 and X x /2 = X z /2 = 15. In 
the initial state, the function q (A) follows a power- 
law with index —2a (see Section 3.1 and Fig. 3). 
In the final equilibrium state, the power-law index 
varies with A, but q(A) is always a monotonically 
decreasing function. As discussed in Section 3.2 and 
shown in Fig. 3, the slope of q in the initial state 
of the axisymmetric model is steeper than in the 
Cartesian model. In fact, qo reaches a maximum 
value of 10 7 ' 3 in the axisymmetric model, but only 
10 6 ' 2 in the Cartesian model. Note that the range 
of the g-axis in Fig. 3 was taken somewhat differ- 
ent in order to facilitate comparison of the shape of 
the curves. The behaviour of q versus $ in the fi- 
nal equilibrium state is similar as q versus A (see 
Fig. 3), although the new geometry introduces some 
remarkable differences in the shape of the function 
q(&). For instance, the change of q, relative to its 
initial value, in the range 0.6 < log 10 < 1.5 is 
smaller than in the same range of A in the Carte- 
sian model. The similar appearance of the functions 
q(A) and q($) should not lead us to the false con- 
clusion that the source terms of the elliptic differen- 
tial equations, namely Q xy (x,y) and Q z r(z, R), are 
also similar. Figure 4 displays Q xy at two different 
cuts along x — and x = A max , and Q z r along 
z = and z = Z max . The fractional change of Q xy 
from its initial value is always larger than the change 
of Q z r- Both Q xy and Q z r may present two local 
maxima along a cut through the maximum heights 
of the magnetic field lines, i.e. along x = A max and 
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Fig. 5. Color map of the magnetic flux <&(z,R) in nat- 
ural logarithmic scale for the final state for A 2 /2 = 12 
(in dimensionless units) and different values of ot. The 
radial gravitational acceleration qr was assumed to be 
constant. 
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Fig. 6. Contour plots of the density (color) and magnetic field lines (solid lines) for the Cartesian model (left) and the 
axisymmetric model (right), under constant gravity, for half the wavelength of 9, 12 and 15 in units of c? a /g (from top 
to bottom) and a = 1. The color bar gives the correlation between the colors and the natural logarithm of the density. 
The number on some curves is the value of A in the Cartesian model or the value of $ in the axisymmetric model. 
To make the comparison of the deformation of the magnetic field lines easier, we plot the same magnetic lines in both 
models. 



THE PARKER INSTABILITY IN FILAMENTS 



11 




Fig. 7. Column density map for the initial (top panel) 
and final axisymmetric state (bottom panel) along a line 
of sight perpendicular to the axis of symmetry (z-axis). 
We used a = 1, \ z /2 = 15 and uniform radial gravity. 
The units of surface density and length are po{Q)(? a /g and 
Cs/g, respectively. Note that the scale is linear. 



z = Z max . The relative amplitude of these maxima 
varies with the wavelength of the perturbation. 

Figure 5 shows contour plots of the magnetic flux 
in the final axisymmetric state for a disturbance with 
X z /2 = 12 and three different a values (a = 0.7, 1.0 
and 1.4). We see that the maps are a scaled version 
to each other by the factor 1 + a. Since it is common 
for astrophysical purposes to assume equipartition 
between magnetic field and gas pressures, we will 
restrict ourselves to the case a — 1 hereafter. 

Figure 6 exhibits contour plots of the density 
(color) and magnetic field lines in the final stable 
state of the axisymmetric model for three different 
horizontal wavelengths (A z /2 = 9,12 and 24). For 
comparison, density and magnetic field lines are also 
shown for the Cartesian model. Although the ver- 
tical size of the computational domain is larger in 
the axisymmetric case, we show the same box in all 
the cases for ease in visualization. When the sys- 
tems are disturbed with the same wavelength, it is 
clear from Fig. 6 that the field lines and the den- 
sity contours become more deformed in the Carte- 
sian model than under axial symmetry. It turns out 



that a cylindrical initial configuration is more rigid to 
perturbations and the magnetic field is less buoyant 
than a plane-parallel layer. Interestingly, the same 
degree of deformation of magnetic field lines found 
for X x /2 — 9 in the Cartesian model, can be achieved 
in the axisymmetric case but for a longer wavelength 
of X z j2 = 12. Likewise, the level of radial buoyancy 
of the magnetic field for X z /2 = 15 is comparable to 
the level of vertical buoyancy for X x /2 = 12 in the 
Cartesian model. 

The column density map of a filament with a = 1 
and whose axis of symmetry lies in the plane of the 
sky is given in Fig. 7 in the initial and final equilib- 
rium states with X z /2 = 15. The scaleheight of the 
column density along a condensation in the direction 
perperdicular to z is about twice the scaleheight in 
the wings of the condensations. As a consequence of 
the less efficient drainage of gas into the magnetic 
valleys, the column density enhancement in the gas 
condensations is rather modest. In fact, the maxi- 
mum enhancement in the column density over a line 
of sight with zero impact parameter, is smaller than 
in the Cartesian model (see Fig. 8). The column den- 
sity may increase up to a factor of 7/4 in Cartesian 
geometry (see also Kim et al. 2004) but only by a 
factor of 5/4 in the axisymmetric model. 

The result that a cylindrical configuration is more 
rigid than a plane-parallel layer has, a posteriori, 
a simple geometrical explanation as follows. Since 
magnetic forces act only perpendicular to the field 
lines, gas pressure gradients must balance the grav- 
itational forces along a field line regardless the ge- 
ometry of the problem. The plasma must slide down 
along the magnetic field lines until it reaches a new 
configuration of pressure equilibrium. The point is 
that, because of the radial convergence of the flow 
streamlines in the axisymmetrical model, we need 
less displacement of gas in the axisymmetric case 
than in the Cartesian model to achieve the same 
pressure gradient, as soon as the gas is isothermal, 
and, henceforth, less deformation of the magnetic 
field lines. 

Since the separation between magnetic field lines 
in azimuthal sections do not reflect the magnetic 
strength in the axisymmetric case, it is illustrative to 
compare the magnetic pressure to make a fair com- 
parison between the results in Cartesian and axial 
geometries (see Fig. 9). Given a certain wavelength 
of the perturbation, the magnetic pressure maps de- 
pend on the geometry. Because of the stronger defor- 
mation in the Cartesian scenario, the magnetic pres- 
sure force and the magnetic stresses along the valleys 
of the magnetic field lines, i.e. at the cut x = 0, are 
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Fig. 8. Column density for the final states as a function 
of the horizontal coordinate for A/2 = 15 and constant 
gravity. In the Cartesion model, the column density is 
S(x) = f^° aa p(x,y)dy. In the axisymmetric model, we 
plot the column density along a line of sight perpendicu- 
lar to the axis of symmetry and with null impact param- 
eter, i.e. = p(z, R)dR. The column density at 
the initial states is drawn as a dotted line. 



larger than along z = 0. The distribution of mag- 
netic pressure in the axial model with A 2 /2 = 12 is 
quite similar to the magnetic pressure configuration 
in the Cartesian case with X x /2 = 9. 

The ratio of the magnetic to gas pressures is an 
indicator of the efficiency with which gas flows along 
field lines. As the increase of gas density in the mag- 
netic field valleys is due primarily to the drainage 
along the field lines (see also M74), we expect that 
if the drainage is less efficient in the axisymmctri- 
cal case, the variation of a in the final state rela- 
tive to the initial values should be also smaller. In 
Fig. 10, we plot the distribution of the magnetic-to- 
gas pressure ratio along x = and x = X max for the 
Cartesian model, and along the same cuts z = and 
z = Z max for the axisymmetric model in the final 
state, for a perturbation with half the wavelength 
of 15. Although the general shape of the curves is 
rather similar in both geometries, the a values in the 
range 10 < y < 22 at x — A max are almost a fac- 
tor 10 larger than in the axisymmetric model at the 
same region (i.e. 10 < R < 22 and z = Z max ). 

4.2. Non-uniform gravitational acceleration 

In a infinite self-gravitating isothermal filament, 
the gravitational acceleration gn at large R decays as 
1/R, whereas it increases linearly at small R. There- 
fore, the assumption that qr is constant is only valid 
at some intermediate range in R. To close our anal- 
ysis, in this subsection we present the final equilibria 
for axially-symmetric configurations in a more real- 
istic external gravitational field. 

According to the results of Kim & Hong (1998) 
for a Cartesian model, the growth time is reduced by 
almost an order of magnitude and the length scale 
of the most unstable mode by factors of 4-8 when a 
more realistic gravity model is used because in such 
a realistic disk the magnetic field strength decreases 
rapidly with vertical distance which promotes the 
development of strong magnetic buoyancy. We ex- 
pect, therefore, that the length and timescales of the 
Parker instability in a filament will be modified if a 
more realistic gravity model is adopted. 

In order to have a model as self-consistent as 
possible, we will assume that the external gravita- 
tional potential is dominated by a collisionless (stel- 
lar or dark matter) component with isothermal one- 
dimensional velocity dispersion o+. The gravita- 
tional potential created by an isothermal filament 
is given by 



^(i?) = 2a> 1 



R- 



(31) 




Fig. 9. Magnetic pressure maps, normalized to its value on the z-axis in the initial state, in the six final states of Figure 
6. 
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where H± is a parameter that specifies the scale 
height of the collisionless component satisfying 
4irGpo : +(0)H 2 = crl, where po,*(0) is the density of 
collisionless matter at R = (e.g., Ostriker 1964). 
As already anticipated, the gravitational accelera- 
tion gn(R) increases linearly with R at R <C y/8H±, 
reaches a maximum at \/8H+ and decreases as 1/R 
at R \/8H*. Under this gravitational field, the 
density and magnetic field in the unperturbed state 
are: 



Po(R) _ B 2 (R) 



= 1 



R 2 



(32) 



Po(0) B 2 (0) ^ ' 

where s is defined by s = 2a 2 /[(l+a)c 2 .]. If we define 
the scale height of the gas, H, as the radial distance 
at which the density of the gas is reduced to 1/e of 
its value at R= 0, then H 2 = 8(exp(l/s) - l)H*. 

To illustrate the effect of using a non-constant 
gravitational acceleration, we useed a = 1 but dif- 
ferent values of s. We ran models where the up- 
per cap was located at the distance i? max at which 
the density decreases to (1 — 5) x 10~ 6 of its value 
at R = 0. We have empirically derived that the 
larger the s- value is, smaller the critical wavelength 
becomes. For s = 4, perturbations with A z /2 = 
3.85iJ* = 2.55H are marginally unstable, whereas 
they become marginally unstable at X z /2 = 1.5H 
when s = 8. In terms of H, this result implies that 
the critical wavelength for s = 4-8 is smaller by al- 
most a factor of 1.5-2.5 than that for the model with 
constant gravity (see Section 4.1) and, thereby, the 
growth timescale is shorter for larger s-values. Fig- 
ure 11 shows the final equilibrium state for s = 4 
and A 2 /2 = 5.3If. This A z -value is close to the hor- 
izontal wavelength with maximum growth rate. By 
comparing with Figure 6, we note that although the 
global appearance is similar to the case with constant 
gravity and X z /2 = 15, the field lines appear slightly 
more deformed. The main effect of a non-constant 
gravitational acceleration strength is to shorten the 
separation of the magnetic valleys. The projected 
surface density that an observer would see if the 
axis of the filament lies in the plane of the sky is 
shown in Fig. 12 for a model with s = 4. The en- 
hancement in the column density over a line of sight 
passing through the core of a condensation is 2.1 in 
this model. This factor is appreciably larger than 
that obtained in axisymmetric models with constant 
gravity but similar to that found in Cartesian models 
with uniform gravity. 

5. DISCUSSION AND CONCLUSIONS 
The Parker instability dictates that longitudinal 
magnetic field lines that give some support against 



gravity are always unstable to undular perturba- 
tions. It has long been realized that the Parker in- 
stability, which is interesting in itself, may play a sig- 
nificant role to understand the formation of massive 
clouds in the interstellar medium, probably aided by 
thermal and gravitational instabilities (Mouschovias 
et al. 2009). The Parker instability can also play a 
role in the evolution of filamentary structure, such as 
filamentary clouds, at least at the early stages of the 
instability. While the properties of the final states 
are well documented for Cartesian models (e.g., M74; 
Basu et al. 1997; Kim et al. 2001), it was unclear how 
they depend on the adopted geometry. In order to 
fill this gap, our primary goal was to investigate the 
nonlinear outcome of such instability in an axisym- 
metric model where the initial equilibrium configu- 
ration consists of a infinite cylinder in the presence 
of a longitudinal magnetic field. We quantified the 
level of buoyancy of the magnetic field and drainage 
of the gas to promote density enhancements, as com- 
pared to the Cartesian models, which are currently 
being used to study the formation of substructure in 
the interstellar medium. We focused on understand- 
ing the physics of the Parker instability in filaments, 
rather than on making detailed comparisons with ob- 
servations. This may be very useful as a first step to 
interpret full numerical simulations. 

In the Cartesian model, the potential vector func- 
tion A obeys a nonlinear elliptic equation, while 
that for the axisymmetric model is the magnetic flux 
$. However, the differential operator as well as the 
source terms are different. In order to gain some 
insight on the different nature of the equations, we 
have compared the form of q(A) and q($) at the ini- 
tial and final states in each model. The source terms 
depend directly on the derivatives of these functions. 
However, in spite that the aspect of these functions 
looks like similar at a glance, the source terms Q xy 
and Q z r turn out to be absolutely different in the 
final state. This is a consequence of the nonlinearity 
of the differential equations. 

We have focused first on models with uniform ex- 
ternal gravity. We found that the axial model is more 
rigid than the Cartesian layer, in the sense that the 
magnetic field is less buoyant and the drainage of gas 
less efficient. In fact, in the axisymmetric model the 
matter shrinks in radius and generates a convergent 
flow that facilitates to achieve pressure equilibrium 
with less deformation of field lines. In order to have 
the same level of deformation of the field lines than 
in the Cartesian model with half the wavelength of 
9c 2 /g, we need a perturbation of half the wavelength 
of 12c 2 / g, for a = 1. Under uniform gravity, we find 
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Fig. 10. Decadic logarithm of the ratio between the mag- 
netic and gas pressures in the final state under constant 
gravity condition. The two solid lines give a values along 
the lines x — and x = A max = 15 for the Carte- 
sian model and the dashed lines give a along z — and 
z — Z max = 15 in the axisymmetric case. 

that a factor 5/4 enhancement in column density can 
be obtained in the axisymmetric model for the more 
unstable wavelength, which is modest as compared 
to the factor of 7/4 found in the Cartesian model. 
At the position of the wings of the condensations, 
the ratio of the magnetic-to-gas pressures becomes 
larger than 1 beyond a certain radius because of the 
evacuation of gas in these regions. However, it is 
still approximately 10 times smaller than in Carte- 
sian geometry. 

We have examined a more realistic choice to de- 
scribe the variation of the external gravitational po- 
tential in axisymmetric models. The main effects of 
a non-uniform gravitational acceleration is to mod- 
ify the separation of magnetic valleys and thereof 
the timescale of the Parker instability. A factor of 2 
enhancement in column density over its initial equi- 
librium value can be found in axisymmetric mod- 
els with non-uniform radial gravity. Hence, unless 
this density enhancement is sufficient to initiate ther- 
mal or gravitational instabilities, the Parker insta- 
bility should be thought of as setting the stage upon 
which other small-scale processes will be modeling 
the structure and evolution of the filament. 
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Fig. 11. Distribution of the volume density and magnetic 
field lines in the axisymmetric final equilibrium state in 
a non-uniform radial gravitational field with s = 4. In 
order to facilitate comparison with Fig. 6 in which H = 2, 
distance is measured in units of H/2. Field lines (solid 
curves) were chosen so that they have the same magnetic 
strength in the initial state than those displayed in Fig. 6. 
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